#Alexander F. Gazmararian
#afg2@princeton.edu

#Load packages
library(readxl)
library(tidyverse)
library(here)
#Load data
g <- read_xlsx(here("data", "input", "fatality", "msha_coalfatalities.xlsx"))
#Prepare data
t1 <- g[, c(1:3)]
t2 <- g[, c(4:6)]
t3 <- g[, c(7:9)]
t4 <- g[, c(10:12)]
t5 <- g[, c(13:15)]
n <- c("year", "miners", "fatalities")
names(t1) <- n
names(t2) <- n
names(t3) <- n
names(t4) <- n
names(t5) <- n
df <- rbind(t1, t2, t3, t4, t5)
p.labor <- df %>%
  filter(year>1972 & year <=2020) %>%
  ggplot(aes(x=year, y = fatalities/miners*10000)) +
  geom_point() +
  stat_smooth() +
  theme_bw(base_size = 14) +
  labs(y="Fatalities per 10,000 workers",x="Year") +
  scale_x_continuous(breaks = seq(1972, 2020, 8)) +
  theme(
    panel.grid = element_blank()
  )
#Figure B1
ggsave(
  p.labor,
  filename = here("output", "figures", "si_fig_B1_minefatalities.png"),
  dpi = 300,
  width = 5.5,
  height = 2.9,
  scale = 1.5
)
